Żeby napisać program komputerowy obliczający takie ciągłe przybliżanie terenu, musimy wykonać następujące czynności.
Po pierwsze, musimy opisać nasz problem w sposób formalny, matematyczny. W szczególności musimy wybrać matematyczną definicję obiektu - bitmapy, oraz matematyczną definicje naszego ciągłego opisu terenu reprezentowanego przez bitmapę.
Rozsądnym rozwiązaniem wydaje się powiedzenie iż bitmapa jest funkcją określoną na obszarze
\( \Omega = [1,maxx]\times[1,maxy] \ni (x,y) \rightarrow BITMAP(x,y) \in [0,255] \)
gdzie poprzez \( \Omega \) oznaczamy cały obszar, na którym rozpięta jest nasza bitmapa.
Z kolei nasza ciągła reprezentacja świata będzie reprezentowana przez funkcję \( u \) o wartościach rzeczywistych.
\( \Omega = [1,maxx]\times[1,maxy] \ni (x,y) \rightarrow u(x,y) \in [0,255] \)
Chcemy by nasza funkcja była "smukła" i ciągła. Matematycznie zapisujemy warunek że nasza funkcja będzie klasy \( C^1 \), oznacza to że w każdym punkcie będzie na tyle smukła żeby dało się policzyć jej pochodne w kierunkach prostopadłych do brzegów naszego obszaru. Oznacza to praktycznie iż w każdym punkcie możemy do wykresu naszej funkcji przyłożyć "linijkę" prostopadle do jednego z brzegów naszego obszaru, i zmierzyć kąt pomiędzy tą linijką a podstawą (powierzchnią płaską rozpostartą na wysokości zerowej). Pochodna to przecież nic innego jak tangens tego kąta. Innymi słowy nasza funkcja nie będzie miała "załamań" ani przerw, na których nie byłoby wiadomo jak przykładać naszą linijkę. Na załamaniach linijka ta przeskakiwała by z jednej pozycji na drugą, a w przypadku dziury w ogóle nie byłoby wiadomo jak ją przyłożyć. Oczywiście możliwość mierzenia pochodnej (przykładania linijki) w dwóch kierunkach prostopadłych do brzegu obszaru oznacza również możliwość przykładania linijki i mierzenia pochodnych kierunkowych w dowolnych innych kierunkach nie prostopadłych do brzegu. Możemy więc gładko przemieszczać się po powierzchni wykresu takiej smukłej funkcji, w dowolnych kierunkach.
Jak uzyskać taką ciągłą smukłą funkcję? Musimy zdecydować jak nasza funkcja zostanie skonstruowana. Na przykład, możemy obszar, na którym rozpięta jest bitmapa podzielić na pewne elementy, i na tych elementach zdefiniować zbiór wielu smukłych funkcji, z których następnie "skleimy" naszą funkcję \( u \).
Wyobraźmy sobie że na wysokości odpowiadającej wysokości zerowej (takiej której odpowiada wartość piksela zero) budujemy płaską dwuwymiarową siatkę, której liczba kwadratowych elementów jest dowolna. Oczka te zgodnie z przyjętą nomenklaturą nazwiemy elementami skończonymi, ponieważ każdy z tych elementów posiada ograniczony skończony obszar. Elementów tych może być mniej niż pikseli w bitmapie, i wtedy nad każdym takim elementem rozpiętych będzie kilka pikseli. Granice pomiędzy naszymi elementami nie muszą zgadzać się z granicami pikseli. Mogą one być dowolnie zdefiniowane na płaskiej powierzchni. Może również być tak, że naszych elementów jest więcej niż pikseli, i wtedy w każdym pikselu znajdować się będzie wiele takich elementów. Załóżmy jednak że naszych elementów jest mniej niż pikseli. Tworzą one regularną siatkę o \( N_x*N_y \) elementach skończonych.
Teraz, na każdym takim elemencie definiujemy smukłą funkcję. Możemy w tym celu posłużyć się funkcjami B-spline.
Funkcje te zostały po raz pierwszy wprowadzone przez amerykańskiego matematyka rumuńskiego pochodzenia Isaaka Jakuba Schoenberga [1].
Funkcje B-spline są powszechnie stosowane w modelowaniu i symulacjach komputerowych dzięki rosnącej popularności dziedziny zwanej analizą izogeometryczną rozpowszechnianej przez prof. T.J.R. Hughesa.[2].
Ideą tych metod jest zastosowanie rodzin funkcji B-spline do obliczeń za pomocą metody elementów skończonych.
Funkcje te oznaczamy \( B_{i,j;2}(x,y) \), gdzie \( i \) oraz \( j \) oznaczają numeracje naszych funkcji, a 2 oznacza iż są to wielomiany kawałkami drugiego stopnia, klasy \( C^1 \).

Rys. 1 przedstawia trzy przykładowe dwuwymiarowe funkcje B-spline rozpięte na dwuwymiarowej siatce. Każda taka dwuwymiarowa funkcja B-spline powstaje poprzez wybranie i przemnożenie przez siebie dwóch jednowymiarowych funkcji B-spline, jednej wybranej ze zbioru jednowymiarowych funkcji B-spline rozpiętych wzdłuż poziomego brzegu siatki, i drugiej wybranej ze zbioru jednowymiarowych funkcji B-spline rozpiętych wzdłuż pionowego brzegu siatki. Zbiory te nazywamy bazami jednowymiarowych funkcji B-spline.
Te jednowymiarowe funkcje B-spline oznaczamy z kolei \( B^x_{i;2}(x) \) oraz \( B^y_{j;2}(y) \) gdzie zmienne \( x \) oraz \( y \) identyfikują kierunki (osie układu współrzędnych) wzdłuż których określone są nasze B-spline'y (oś pozioma \( x \) oraz oś pionowa \( y \)), \( i \) oraz \( j \) oznaczają numeracje tych funkcji (którą kolejną jednowymiarową funkcję B-spline wybieramy z takiej jednowymiarowej bazy), oraz 2 oznacza ponownie że są to wielomiany kawałkami drugiego stopnia, klasy \( C^1 \) (czyli że umiemy z nich liczyć pierwsze pochodne).
Następnie wybrane funkcje jednowymiarowe są przemnażane przez siebie, co daje nam smukłą dwumymiarową funkcję B-spline. Taką metodę tworzenia dwuwymiarowych funkcji przez przemnażanie stosownych jednowymiarowych funkcji nazywa się iloczynem tensorowym funkcji.
\( B_{i,j;2}(x,y)=B^x_{i;2}(x)B^y_{j;2}(y) \)
Jest to zilustrowane na Rysunku. Uzyskane w ten sposób dwuwymiarowe funkcje B-spline mają kształt smukłego "pagórka", określonego na dziewięciu sąsiadujących elementach. Najwyższy punkt takiej funkcji - pagórka znajduje się w centrum środkowego elementu. Funkcje te schodzą gładko do wartości zerowej, przyjmowanej na brzegach kwadratu zdefiniowanego przez dziewięć elementów, na których funkcja ta jest określona. Funkcje te zgodnie z przyjętą nomenklaturą nazwiemy funkcjami bazowymi.

Naszą ciągłą aproksymacje terenu uzyskamy w ten sposób, że zsumujemy ze sobą wiele takich smukłych pagórków - B-spline'ów. Każdy z nich zostanie przeskalowany (podniesiony do góry lub dół) tak by w sumie otrzymać ciągłą aproksymację terenu. Jeśli poprawnie dobierzemy wysokości poszczególnych pagórków, dostaniemy wówczas smukłe przybliżenia naszego terenu tak jak przedstawiono to na Przykładowy problem dwuwymiarowej projekcji bitmapy-Rys. 2.
Powstaje teraz pytanie, w jaki sposób dobrać wysokości, do których wyciągniemy nasze funkcje bazowe. Pierwsza metoda, która przychodzi nam do głowy to wybrać wartość piksela z \( BITMAP(x,y) \) znajdującego się dokładnie w najwyższym punkcie B-spline (na środku pagórka). Niestety metoda ta ma kilka wad. Po pierwsze, jeśli naszych funkcji bazowych B-spline jest mniej niż pikseli, to wówczas ignorujemy wszystkie sąsiednie piksele znajdujące się w obszarze naszego B-spline'a, wybierając jedynie jedną wartość ze środka obszaru. Możliwe że nasza bitmapa posiada pewne zaburzenia i że trafimy akurat na lokalny odskok będący błędem pomiaru, lub że trafimy akurat w lokalną dziurę w ukształtowaniu terenu, lub w lokalne drzewo lub budynek który akurat zakłócił pomiar ukształtowania terenu. Po drugie, zauważmy że nasze funkcje bazowe rozpościerają się na kwadracie dziewięciu elementów. Skoro każdą taką funkcję B-spline skojarzyć można ze środkiem swojego elementu, oraz każda z tych funkcji rozpościera się na dziewięć sąsiednich elementów, oznacza to iż na każdym elemencie rozpostartych jest w sumie dziewięć takich funkcji, oraz że sąsiadujące funkcje zachodzą na siebie. Jeśli więc rozciągalibyśmy funkcje tak aby ich punkt maksymalny pokrywał się z centralnym pikselem, oraz zsumowalibyśmy te wszystkie funkcje razem w celu uzyskania naszej globalnej funkcji \( u \), to wówczas na każdym elemencie, nawet w punkcie centralnym, nasza wynikowa aproksymacja (suma tych funkcji) byłaby wyżej niż nasz centralny piksel, dlatego iż sąsiednie osiem funkcji również byłyby niezerowe na danym elemencie i podnosiłyby one wartość naszej aproksymacji w tym punkcie do góry.